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1 Introduction 

At a second-order phase transition different types of order of a physical system meet 
one another, giving rise to universal critical properties which are independent of the 
microscopic details of the system. This kind of universality is an extremely successful 
concept in characterising equilibrium states of matter and classifying different phe- 
nomena in terms of just a few classes governed by the same critical properties. The 
appearance of an ordered state which does not possess the full symmetry deriving from 
the conservation laws obeyed by the system's dynamics is called spontaneous symme- 
try breaking. The order induced by this symmetry breaking has long been known to 
allow for the appearance of (quasi-)topological defects such as solitons or vortices, 
objects we are particularly interested in here. 

In these notes we leave equilibrium systems aside and aim at sketching a pic- 
ture of related concepts far from thermal equilibrium. We consider the example of 
nearly coherent Bose gases brought far out of equilibrium and discuss their behaviour 
in view of connections between universal properties, (quasi-)topological field config- 
urations and turbulent dynamics. We demonstrate that the isolated Bose gas, on its 
way back to thermal equilibrium, can approach metastable non-equilibrium configu- 
rations and spend a long time in their vicinity. In such configurations, which have 
been termed non-thermal fixed points, the system shows universal long-range proper- 
ties manifest through scaling, i.e., self-similar correlations. The time evolution near 
such fixed points is demonstrated to undergo critical slowing down. The spatial field 
pattern, at the same time, is characterized by the appearance of defects and domain 
formation whose geometry gives rise to the particular scaling laws seen in the cor- 
relation functions. We obtain an overall picture which connects well-known concepts 

1 Prepared for the proceedings of the Summer school: Strongly interacting quantum systems out 
of equilibrium, held 30 July-24 August 2012 at Ecolc dc Physique des Houchcs, Lcs Houches, France, 
to which T.G. contributed a seminar talk. 
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for describing universal dynamics such as wave-turbulence, supcrfluid turbulence, and 
(quasi-)topological excitations. This allows to bring together an excitingly wide range 
of concepts and methods with an excitingly wide spectrum of applications. Beyond the 
immediate implications for simple low-energy degenerate quantum gases, phenomena 
such as topological configurations in solids, in soft matter, the dynamics of the quark- 
gluon plasma created in heavy-ion collisions, or the reheating of the post-inflationary 
universe come in sight. Vice versa, an ultracold quantum gas offers itself as a kind 
of the often discussed 'quantum simulator' for universal dynamics of systems which 
are more difficult to access experimentally, like some of the above. Not the least, the 
extension to fcrmionic and gauge fields should bear many new interesting aspects. 

Time-evolution far from equilibrium of a system with many degrees of freedom 
is often characterized by the appearance of widely different scales. Usually, a period 
of fast motion is followed by a much longer period of much slower motion. In the 
same way a separation of spatial scales allows descriptions in terms of statistical con- 
cepts such as hydrodynamics or renormalization-group theory. Taking into account 
the well-developed concepts and concrete results for equilibrium and near-equilibrium 
systems (cf., e.g., Hohenberg and Halperin, 1977), it has to be expected that far- 
from-equilibrium time evolution will show universal behaviour, with fixed points or 
partial fixed points being just a small subclass of the possible phenomena. Recently, 
the discussion of the possible things waiting to be discovered has increased in intensity, 
in particular in the context of possible 'prethermalisation' phenomena (Berges et al., 
2004; Bonini and Wetterich, 1999; Aarts et al., 2000; Gasenzer et al, 2005; Berges and 
Gasenzer, 2007; Barnett et al, 2011; Kitagawa et al, 2011; Gring et al., 2012; Kollar 
et al., 2011). Such dynamics could be considered to hint to a wider class of universal 
phenomena which will allow to classify and more deeply understand non-equilibrium 
physics. 

2 Strong wave turbulence and vortical flow 

In analogy to equilibrium phase transitions and criticality in driven systems it has been 
proposed that transient stationarity arising in the time evolution of an initially strongly 
perturbed quantum many-body system can reflect the existence of non-thermal fixed 
points (Berges et al., 2008). At the same time, these fixed points were numerically 
demonstrated to occur in the field evolution imposed by an 0(A)-symmetric, relativis- 
tic, non-linear scalar model: By means of parametrically resonant oscillations of the 
field expectation value, and non-linear amplification, a broad range of modes can be ex- 
cited initially. While falling back to equilibrium their occupation number spectra shows 
scaling behaviour. Remarkably, the respective power-law exponents confirmed analyt- 
ical predictions in the infrared domain of long-wavelength excitations where standard 
descriptions in terms of Boltzmann-type kinetic equations break down. Found when 
searching for stationary scaling solutions of non-perturbative dynamic Dyson equa- 
tions, with the help of Zakharov integral transformations known in wave-turbulence 
theory (Zakharov et al., 1992; Nazarenko, 2011), these infrared power laws are inter- 
preted to constitute the previously mysterious strong wave turbulence which had been 
considered out of reach of kinetic theory (Scheppach et al., 2010). 
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In the following we identify this strong wave turbulence, considering the example of 
a non-relativistic, superfluid Bose gas in two spatial dimensions, with the appearance 
of quantised vortices. On scales considerably smaller than the mean distance between 
vortices, the velocity field v(r) associated with the rotational flow decays as 1/r with 
growing distance r from the nearest vortex core. As a result, the angle-averaged ki- 
netic energy distribution ~ mv 2 /2 of Bosons gives rise to the single-particle momen- 
tum spectrum n(k) <~ k~ A identical to that predicted for strong wave turbulence by 
Schcppach et al. (2010). Hence, the non-thermal fixed point appears to correspond to 
a configuration bearing a dilute ensemble of vortices. In this section we give evidence 
for this interpretation which allows to conjecture a deep link between the extended 
kinetic-theory picture of wave turbulence on the one side and the theory of non-linear 
(quasi-)topological field configurations and superfluid turbulence on the other. 

2.1 Non-thermal fixed points and weak wave turbulence 

Most generally, a non-thermal fixed point can be defined as a metastable state of a 
many-body system (Bonini and Wetterich, 1999; Aarts et al., 2000; Berges et al., 2008). 
While one may include open systems into the discussion, e.g. driven ones and those 
experiencing dissipation (Diehl et al., 2008), we will, in the following, restrict ourselves 
to closed systems. In analogy to fixed points in scaling flows, one furthermore consid- 
ers in particular solutions with power-law behaviour of correlation functions (Berges 
et al, 2008; Berges and Hoffmeister, 2009; Scheppach et al, 2010). While, precisely 
at the fixed point, the system is metastable and characterised by scaling correla- 
tion functions in the infrared limit of infinitcsimally slow modes, it will vary in time 
away from the fixed point and show power-law correlations within a finite scaling 
regime. From the physics of turbulence, it is well known that such states occur as a 
consequence of local conservation laws in momentum space. Such unintuitive locally 
conserved currents in momentum space instead of position space imply the possibil- 
ity of gain- and loss-less transport processes between different scales, giving rise to 
so-called cascades. The most prominent example of this kind is fully developed classi- 
cal fluid turbulence. It comprises a quasi-stationary flow of kinetic energy from large 
to small spatial scales, i.e., from low to high momenta (Richardson, 1920). The en- 
ergy is fed in, e.g., by a stirrer, at a large scale and finally dissipated into heat at 
the microscopic scale defined by the fluid's viscosity. The corresponding energy spec- 
trum exhibits the famous Kolmogorov-Obukhov five-third scaling of the radial energy 
distribution E(k) ~ fc~ 5 / 3 (Kolmogorov, 1941; Obukhov, 1941). 

Considering a dilute Bose gas one has to take into account its compressibility, 
allowing for collective sound excitations of the particles. In addition to the density 
of an incompressible fluid the system is thus characterised by the dispersion relation 
between momentum and energy of its excitations. This allows for the techniques of 
wave turbulence to be invoked when looking into turbulence phenomena of a Bose gas 
(Zakharov et al, 1992; Nazarenko, 2011). The mathematically best-controlled case is 
that of weak wave turbulence which rests on the analysis of stationary solutions of 
Boltzmann-type kinetic equations. Within a certain range of momenta k and times t, 
for not too strongly excited systems, the quantum Boltzmann equation (QBE) 
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d t n(k,t) = J(k,i), (1) 

J(k,t) = y d d pd d gd'V|T kpqr | 2 ( 5(k + p-q-r) ( 5K + Wp -w q -w P ) 

x [(n k + l)(n p + l)n q n r - n k n p (n q + l)(n r + 1)], (2) 

well describes the time-evolution of the momentum-mode occupation numbers n k = 
n(k, t) = ($t(k, i)$(k, i)) of an interacting degenerate Bose gas. Here, $ denotes the 
quantum field operator describing the Bose system, [3>(x, t), <f?' (y, t)] = 5(x — y), while 
all other equal-time commutators vanish. We consider only two-to-two elastic collisions 
quantified by the T-matrix elements which for dilute, weakly interacting atomic gases 
reduce to a single quantity, the s-wave scattering length a, i.e., |T kpqr | = g = const, x a. 

Zeroes of the scattering integral I(k) correspond to fixed points of the time evo- 
lution within the regime of applicability of the QBE (Zakharov et ai, 1992). Most 
prominent amongst these are the thermal fixed point corresponding to the system in 
thermal equilibrium, ra k = {cxp[/3w(k)] — and the trivial fixed point n k = const. 

At both fixed points the scattering integral vanishes and n(k, t) becomes independent 
of t. Note that both, the trivial and the Bosc-Einstein distribution (in the Rayleigh- 
Jcans regime), taking w(k) <~ k 2 , show a power-law behaviour of the form n k ~ k^^ 
with C = and ( = 2, respectively. 

The theory of weak wave turbulence (Zakharov et ai, 1992) allows to analytically 
derive further, non-thermal fixed points at which the occupation number n k obeys a 
scaling law of the form n k <~ k~^ with, in general, ( ^ 2. As in classical turbulence one 
expects that universal scaling appears within a certain regime of momenta, the inertial 
range. According to this picture, outside the scaling regime excitation quanta enter the 
system from an external or internal source and/or leave it into a sink, whereas there are 
no sources and sinks within the inertial interval where the quanta are transported from 
momentum shell to momentum shell. This process is described by a continuity equa- 
tion in momentum space, with a momentum-independent, radially oriented current 
vector. 2 A central aspect of weak-wave-turbulence theory is that the quantum Boltz- 
mann equation can be cast into different such equations (Zakharov et ai, 1992), for the 
radial densities of particles, N(k) = (2k) d ^ 1 Trn(k), and energy, E(k) = (2k) d ~ 1 ire(k), 
e(k) — Lo(k)n(k), 

d t N(k,t) = -d k Q(k), (3) 

d t E(k,t) = -d k p(k). (4) 

Taking either the radial particle current Q(k) — (2k) d ~ 1 nQk(k) or the energy current 
P(k) = {2k) d - 1 irP k {k) to be independent of k, one derives different scaling exponents. 
The resulting exponents 3 arc 



2 Note that justification of this assumption, i.e., locality of the transport, needs to be checked for 
each particular wave-turbulent solution, cf., e.g., Zakharov et al. (1992). 

3 The superscript UV (ultraviolet) in Eq. (5) refers to the regime of large momenta where the 
description in terms of a kinetic equation is expected to be accurate. 
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Fig. 1 2PI diagrams of the loop expansion of T2 [G] . (a) The two lowest-order diagrams of 
the loop expansion which lead to the quantum Boltzmann equation. Black dots represent the 
bare vertex ~ g8(x ~ y), solid lines the full propagator G(x, y). (b) Diagram representing the 
resummation approximation which, in the IR, replaces the diagrams in (a) and gives rise to 
the scaling of the T-matrix in the IR regime, (c) The wiggly line is the two-point resummed 
vertex function which is represented as a sum of bubble-chain diagrams. 



C^=d-2/3, Cr=d- (5) 

These exponents can be obtained by simple power counting: Combining Eqs. (1) and 
(3) gives the radial relation dkQ(k) ~ k d ~ 1 I(k) which implies that stationarity requires 
k d I(k) to become ^-independent, i.e., scale as k°. Counting all powers of k in I(k), 
Eq. (2), in the wave-kinetic regime where the terms of third order in the occupation 
numbers dominate the scattering integral, this requires n(k) ~ fc~ d+2 / 3 . Analogously 
one infers the exponent (p Y from the balance equation (4) for the energy density 
e(k) ~ k 2 n(k). Despite this simple procedure, the existence of the respective scaling 
solutions has to and can be derived rigorously from the quantum Boltzmann equation 
by means of Zakharov conformal integral transforms (Zakharov et ai, 1992). 

As we will illustrate for our case in Sect. 3.2, the energy flux generically constitutes 
a direct cascade to larger k, whereas the particle flux corresponds to an inverse cascade. 
While the character of the fluxes is entirely determined by the properties of the system, 
it turns out that there must be at least two sinks, where the particles and the energy 
can flow to (Zakharov et ai, 1992). Let us assume that an external source introduces 
particles with energy lo s at some scale k s and that energy and particles leave at scales 
k = and k — K. Calling T s , Tq, and Tk the injection/ejection rates of particle 
number at the respective scales, number and energy conservation imply T s = T + Yr 
and lo s T s = lu T + u>i<T k> respectively. Inverting these conditions, 

LJ K - L0 S LJ S - UJQ 

r = r s , v K = r s , (6) 

shows that for loq -C uj s <^ ujk the particles are ejected at loq whereas energy is 
dissipated at u>k (Gurarie, 1995). 



2.2 Infrared scaling as strong wave turbulence 

Given a positive scaling exponent £ momentum occupation numbers n(k) ~ k~^ grow 
large in the IR regime of small k. Keeping the coupling g fixed, the QBE fails for 
rik > g^ 1 where perturbative contributions to the scattering integral I(k) of order 
higher than g 2 are no longer negligible. To find scaling solutions in the IR, an ap- 
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proach beyond this perturbative approximation is required 4 . This is available through 
quantum-field dynamic equations derived from the two-particle irreducible (2PI) effec- 
tive action or ^-functional (Luttinger and Ward, 1960; Baym, 1962; Cornwall et al., 
1974) which can be expanded in terms of 2PI closed loop diagrams the lowest-order 
ones are sketched in Fig. la. The solid (blue) lines denote the time-ordered Green's 
function G(x, y) — (T& (x)$(y)) c which, in turn, is a solution of the real-time Dyson 
equation. The Dyson equation is derived from the action by use of Hamilton's vari- 
ational principle. It contains a time-evolution equation (1) for the momentum-mode 
occupation numbers n^. As before one considers zeros of the resulting scattering inte- 
gral which can be expressed, within the 2PI approach, in terms of the self-consistently 
determined Green's function G, connected by the bare scattering vertices of the the- 
ory. The scattering integral (2) of the QBE (1) is recovered within the expansion of 
the action up to the two-loop diagrams in Fig. la. 

To describe the IR kinetics one needs to go beyond this approximation. Resumming 
an infinite set of loop diagrams contributing to the 2PI effective action (Bcrgcs, 2002; 
Aarts et al, 2002), also, e.g., (Berges, 2005; Gasenzer et al., 2005; Gasenzer, 2009) 
leads to a non-perturbative, effectively renormalised coupling g e s(k) in the dynamic 
equations (Berges et al., 2008; Berges and Hoffmeister, 2009; Scheppach et at, 2010), 
see Fig. lb and c. In particular, this coupling becomes suppressed in the IR to below its 
bare value g which, in effect, leads to an even steeper rise of the particle spectrum rik 
(Berges et al., 2008). The IR scaling exponents for the radial particle and energy flows 
of a Bose gas in d dimensions which constitute the resulting strong wave turbulence 
were derived by Scheppach et al. (2010) to be 

C^=d + 2, Cp R -d + 2 + z, (7) 

where z is the dynamical scaling exponent accounting for the scaling of the dispersion 
oj(sk) — s z uj(k). From the point of view of its scaling, the T-matrix in the scattering 
integral of the QBE can be replaced by an effective many-body T-matrix, = 
^k+ P q+r- This effective T-matrix is more complex than but scales like 

|7k ff | = ~ \9Ck z - 2 /[l + C'gk d -W]\, (8) 

k = |k|, where C is some constant which fine-tunes the position of the transition from 
UV to IR scaling. The second term in the denominator can be related to the validity 
criterion of the kinetic equation (1) (cf., e.g., Svistunov, 1991) in d dimensions, 

gj\«k>nV)« 1 ^. (9) 

For a scaling distribution n(k) ~ this translates into 2£ldmgk d ~ 2 n(k) <C 1, with 
fid the surface of a unit sphere in d dimensions. For small rik and z = 2 one recovers the 
weak- wave-turbulence case discussed in the previous section. For large n^, the second 



4 For more details we recommend to consult Scheppach et al. (2010) for details of the procedure 
summarised in the following. 
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Fig. 2 Sketch of the single-particle mode occupation number n(k) of a Bose gas in d dimen- 
sions, as a function of the radial momentum k of the system at the non-thermal fixed point, 
as predicted by Scheppach et al. (2010). At large momenta, weak wave turbulence scaling 
is recovered, with the Kolmogorov-Zakharov spectra £ = d — 2/3 (particle flux) and £ = d 
(energy flux). At low momenta, new scaling solutions are predicted, with C, = d + 2 (particle 
flux) and £ = d + 2 + z (energy flux), termed strong turbulence. Colours indicate the expected 
character of an IR particle flux and a UV energy flux. 

term in the denominator dominates which implies a power-law behaviour |T£ ff | 2 ~ 
as a consequence, the modified scaling (7) of rik in the infrared regime 
of small wave numbers. Moreover, the coupling becomes universal in the sense that it 
is now independent of g which is cancelled out by the leading denominator term in 



In Fig. 2, we summarise the non-thermal fixed point scaling predicted with the 
loop-resummed 2PI effective action, for a dilute Bose gas in d dimensions. Direct 
confirmation of the scaling Eq. (7) by integration of the 2PI dynamic equation is 
complicated by the required computational effort. However, as we will show in the 
next chapter, this challenge can be met by use of semi-classical simulations of the field 
equations of motion. 

Before we proceed, let us comment on the possibility to understand the turbu- 
lent scaling taking a renormalisation-group viewpoint (Berges and Hoffmeister, 2009; 
Berges and Mesterhazy, 2012; Philipp, 2012). Turbulence has served, since the semi- 
nal work of Kolmogorov (Kolmogorov, 1941; Obukhov, 1941; Frisch, 1995), as one of 
the first phenomena to develop renormalisation-group techniques out-of-equilibrium. 
The effectively local transport processes in momentum space, which are at the ba- 
sis of turbulent cascades, immediately suggest themselves for a renormalisation group 
analysis. Building on functional renormalisation-group techniques (Wetterich, 1993; 
Berges et al., 2002; Gies, 2006; Pawlowski, 2007) specifically out of equilibrium (Canet 
et al., 2004; Kehrein, 2005; Mitra et al., 2006; Zanella and Calzetta, 2006; Canet and 



(8). 
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Chate, 2007; Gezzi et at, 2007; Jakobs et al, 2007; Korb et al, 2007; Matarrese and 
Pietroni, 2007; Karrasch et al, 2008; Jakobs et al., 2010; Schoeller, 2009; Gasenzer and 
Pawlowski, 2008; Berges and Hoffmeister, 2009; Gasenzer et al, 2010), more refined 
scaling analyses are being developed (Canet et at, 2010; Mathey and Polkovnikov, 
2010). The scaling exponents given in Eq. (7) result as canonical exponents, from power 
counting of flow equations for G and higher-order vertices. Additional flow equations 
which, in particular, account for the interaction effects in the spectral functions, are 
expected to fix potentially relevant anomalous dimensions. We remark that the vor- 
tex picture developed in the following leads to the possibility to conjecture that such 
anomalous scaling will become important in driven systems where interaction effects 
between finitely spaced defects become relevant. 



2.3 Vortices in a superfluid 

Near degeneracy, the strongly occupied low-momentum modes of an ideal Bose gas, 
rak 3> 1, constitute the Rayleigh- Jeans distribution riRj(fc) = 2mT/k 2 of classical 
waves. 5 One can show that also out of equilibrium, strongly occupied modes can be 
described by the dynamics of classical waves (Brewczyk et al., 2007; Berges and Gasen- 
zer, 2007; Blakic et al., 2008; Polkovnikov, 2010). The state of the gas can be defined in 
terms of the Wigner quasi-probability distribution W[<f>, <j)*} for the complex field </>(x) 
and its conjugate momentum <jf (x) at each point in space. In the considered wave- 
classical limit, W is positive definite. Since quantum effects arising from coupling 
to sparsely occupied modes are small, the dynamics of the Wigner function follows a 
classical Liouville equation. In this semi-classical limit, the so-called truncated Wigner 
approximation allows to follow the evolution exactly, within computational errors, by 
evaluating many trajectories evolving according to the classical field equation 



i<9 t 0(x, t) = 



2 



?- + V(x)+g\cf ) ( X ,t)f 



0(x,t). (10) 



While (10) is equivalent in form to the Gross-Pitaevskii equation (GPE) for the quan- 
tum field expectation value ($) the statistical sampling procedure leads to a quasi- 
exact result for the full many-body evolution. Correlation functions are obtained by 
averaging over many trajectories. The non-linear classical field equation (10) has some 
interesting properties relevant for our purposes: It can be mapped to an Euler-type 
hydrodynamic equation, implying the interpretation of the gas dynamics in terms of 
(superfluid) flow. Among the possible solutions of this equation those resembling eddy 
flow and shock waves lead to topologically nontrivial configurations. These include 
vortices, solitons, and related nonlinear stationary states if more than two field com- 
ponents couple to each other, see, e.g. Pismen (1999), Pitaevskii and Stringari (2003). 

Hydrodynamic representation. The polar representation <j) = y/n cxp{if} allows 
to express the particle current j = i(cj>* V0 — 0V<^*)/(2m) = nv in terms of the velocity 
field v = m~ 1 \7tp and the particle density n = |</>| 2 . With this, the GPE (10) can be 



5 In the following we use units where h = kg = 1. 
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rewritten into the continuity and an effective Euler equation for a compressible inviscid 
(i.e. non-viscous) fluid with modified pressure m _1 /i, 

d t n + V • (nv) = 0, (11) 
d t v + (v • V)v = -m _1 V/i, fL = gn- — ^ A\/n. (12) 

As the velocity is a potential field it is irrotational wherever the density is non- 
vanishing and therefore the phase ip of the field <j) wen defined. Small excitations 
of <f) which give rise to small density fluctuations and small velocities are described 
by the linearised version of the hydrodynamic equations (12). At wave lengths larger 
than the healing length £ = prang] -1 / 2 these are collective sound-wave excitations. 

Solitons are quasi-topological one-dimensional solutions of Eq. (10) which travel with 
a fixed velocity but are non-dispersive, i.e., stationary in shape, see, e.g., Zakharov 
and Shabat (1972), Kevrekidis et al. (2008). For positive coupling constant g > 0, the 
solitons are 'dark', i.e., characterised by an exponentially localised density depression 
in the surrounding bulk matter. This and the corresponding phase shift in complex 
field are given by 



K{x,t) = \fn 



(13) 



where x s (t) = xo + ut is the position of the soliton at time t. Depending on the depth of 
the density depression, the dark soliton is either called grey or, for maximum depres- 
sion, black. 7 = 1/vT — v 2 is the 'Lorentz factor' corresponding to the velocity v of the 
grey soliton in units of the speed of sound, v = v/c s — \<f>„(vt, t)\/y/n. Being related 
to the density minimum, v measures the 'greyness' of the soliton, ranging between 
(black soliton, \<j> v (vt,t)\ — v^fn = 0) and 1 (no soliton, \(j) v {vt,t)\ = y/n). Due to the 
interaction with sound, solitons can continuously vanish which means that they are 
not topologically stable. In d > 1 dimensions, solitons decay into vortices (Anderson 
et at, 2001; Brand and Reinhardt, 2002). The energy to create a soliton on top of a 

uniform background is E s = (l — ^ 2 ) 3 ^ 2 4nc s /3 (Pitaevskii and Stringari, 2003). For 
small velocities v <C 1, this gives E s ~ 4nc s /3 — 2nc s v 2 reminiscent of a classical point 
particle with negative mass — 4n/c s . The energy of a soliton monotonously decreases 
with increasing velocity which hints at a dynamical instability. 



Vortices are topologically stable solutions of Eq. (10) in d > 1 dimensions which form 
the superfluid analogies of eddy flows in classical fluids. We recall that the ground- 
state manifold given by the minimum of the effective potential V{4>) = mH 2 + 
of the classical field </>(x) = \J n(x) exp{«<^(x)} requires constant density n(x) = n, 
but is degenerate in the phase </?(x). The true ground state has a constant phase, 
and is therefore called topologically trivial. On the other hand, we can consider field 
configurations which have constant density on the boundary of, e.g., a two dimensional 
volume, but a varying phase. If we use this freedom to evolve the phase angle <^( x ) 
from to 2tt when going around the boundary, we arrive at a topologically nontrivial 
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state. Configurations are topologically distinct, because one can not define a continuous 
function that transforms one into the other. As a consequence of the phase winding the 
phase can not be well-defined at some point inside the volume, and hence, the density 
has to go to zero at that point. The stationary state of Eq. (10) which exhibits these 
properties is called a vortex (Pitaevskii, 1961). Following the phase angle along a closed 
path around the vortex core it continuously varies from up to ip = 27tk, where the 
integer k is called the winding number or circulation. Only a singly quantised vortex 
with k = ±1 is stable. It is described, in polar coordinates centred at the vortex core, 
by the field <p(r, ip) — f(r) cxp{in(p} where f(r) can be chosen real and approaches 
the square root of the bulk density n for large distances r from the vortex core. At 
r — > 0, f(r) <~ r rises linearly. <fi(r, ip) is a stationary solution of Eq. (10), evolving as 
4>(r,t) = <p(r, 0) cxp{— i/it} with /i = gn. 

We remark that the irrotational nature of the velocity field defined in the hydro- 
dynamic formulation of the GPE is restricted to those points where the density is 
non-vanishing. The velocity field of the vortex is, in the polar coordinates used before, 
v(r, ip ) = Ke v /(mr). Its curl is concentrated locally at the vortex core. The compact 
phase p € [0, 2tt), becomes a non-compact velocity potential. As a consequence, vor- 
tex creation or annihilation is not described by the effective Euler equation. In fact, 
due to the Thomson circulation theorem vorticity is locally conserved in an inviscid 
flow (Lesicur, 2008). 

The energy associated to the vortex is extremely non-local: For a singly quantised 
vortex in a 2D homogeneous gas the energy within the volume V — irR 2 grows logarith- 
mically with the radius, E v — 7mra _1 ln (1.46i?/£), see, e.g., Pitaevskii and Stringari 
(2003). In d = 3 dimensions, point vortices extend to vortex lines around which the 
fluid rotates (Pitaevskii, 1961; Pismen, 1999). In the simplest case, a vortex line of 
length L goes straight from one end of the volume to the other. This requires the ex- 
citation energy LE V . Vortex lines can not end inside the medium, but can form closed 
loops of all shapes, e.g. rings, ellipses, also knots. The GPE (10) moreover supports lin- 
ear wave excitations of the position of the vortex lines, so called Kelvin waves (Sonin, 
1987; Kivotides et ai, 2001; Vinen et ai, 2003; Krstulovic, 2012). Higher dimensional 
vortices exist whereby the dimensionality of vortex-core geometry is always d — 2, for 
example leading to vortex surfaces in four dimensions. 

Decomposition of the flow field. A non-equilibrium flow features the presence of 
multiple types of excitations. In order to distinguish longitudinal excitations (sound 
waves) from rotational excitations (vortices), we close this section by discussing a 
decomposition of the kinetic energy density proposed by Nore et al. (1997). The total 
kinetic energy E^ = m J d d x (\ V0(x, t)\ 2 }/2 can be split, E^ n = E v + E q , into a 
'classical' part E v — mf d d x (|y / nv| 2 )/2 and a 'quantum-pressure' component E q = 
J d d x (IVv^l 2 ) / (2m). The radial energy spectra for these fractions involve the Fourier 
transform of the generalised velocities w v = yfnv and w q = V^/n/m, 




(14) 
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which we cast further into occupation numbers ns(k) = k~ 2 Es(k),5 = v,q. Since 
the superfluid velocity v = V<^ is a potential field it does not reveal a transversal 
flow component, V x v = (outside vortex cores). On the contrary, w v is not a 
potential field and the divergence of v at r — > is regularised by the vanishing yjn(r). 
Following Nore et al. (1997) the regularised velocity w v can be further decomposed into 
'incompressible' (divergence-free) and 'compressible' (solenoidal) parts, w v = w;+w c , 
with V • Wj = 0, V x w c = 0, to distinguish vortical superfluid and rotationless motion 
of the fluid. By construction, the generalised velocity w v = /(r)v(r, ip) of a vortex has 
only an incompressible component, since 

V • w v = V/(r) • v(r, <p) + /(r)V • v(r, <p) = . (15) 

The first term vanishes due to the transversal nature of the vortex velocity field, the 
second one equals zero because the superfluid velocity is a potential field. The density 
of incompressible energy |wi| of a vortex is constant up to about one healing length 
distance from the core and then falls off as 1/r. Sound waves are purely compressible 
excitations. In our simulations, their oscillating density and phase profiles will be 
visible as maxima and minima in the compressible energy density in position space. 

2.4 Vortex statistics 

To understand the implications of vortex defects appearing in the dynamical evolution 
of degenerate Bose gases, in particular to make contact to the observables studied in 
the context of wave turbulence we turn to a statistical viewpoint. The point vortex 
model studied in the following was introduced by Onsager (1949). It describes the com- 
plex flow pattern in terms of the statistics of classical point objects with Coulomb-type 
interactions. However, due to the absence of a kinetic energy term in the Hamiltonian 
there is no kinematic transfer of potential into motional energy. The model is con- 
structed as a discrete-vorticity approximation of classical fluid turbulence, but it is 
even more suitable to describe superfluid turbulence consisting of quantised vortices. 

We restrict our discussion to the example of d = 2 dimensions. We have seen 
that an isolated, singly quantised vortex is described by the complex field <p(r, ip) = 
y/n(r)e lip . As the r-dependence of the density n(r) only becomes important at small 
scales on the order of the healing length £ = (2m<?n)~ 1//2 at which in practice thermal 
excitations dominate, we assume n to be uniform. A set of M vortices can be described 
by <j>(x) — Ilf f (/)i(x), where <fo(x) = </>(x — Xj) is the single-vortex field centred around 
Xj. Let us derive the corresponding bosonic single-particle spectrum by considering 
the velocity field v = Vtp/m. We can express the mean classical kinetic energy density 
of the velocity field v(x) of a single vortex as 

£ v (x) = |<|v(x)| 2 )^<| J dV v(x - x') p(x') | 2 ) , (16) 

where p(x) = J2fL\ Ki<5(x — Xj) defines the spatial distribution of vortices with winding 
number Ki = ±1. Here and in the following, (•) denotes an ensemble average over 
different realisations of the classical field 0(x). We derive the low-Zc scaling of n(k) from 
the kinetic-energy spectrum E v (k), the angle-averaged Fourier transform of E v (x), 
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taking into account that at low k, the velocity field v dominates the dynamics, 

n(k) ~ 2mk- 2 E v (k). (17) 

One has, from Eq. (16), 

S v (k)^(|v(k)| 2 )^(|p(k)| 2 |v(k)| 2 ), (18) 

with 

M 

\ P (k)\ 2 = Y,^ eM ^ } - ( 19 ) 



Below the healing length scale k$ — 2tt/£,, the modulus of the velocity field of a single 
vortex scales as |v| <~ k~ x and is radially symmetric. 

To distinguish contributions from vortex-vortex and vortex-antivortex correlations 
we write the distribution p(x) = p v (x) — p A (x) as the sum of distributions /? v (x) = 
Sj=i ^( x — X Y) 01 M vortices and p A (x) = JZ i=1 S(x — x A ) of M antivortices. Hence, 

(\p(k)\ 2 )= J d 2 a;d 2 a;'e jk(x - x ' ) C(x,x'), (20) 

with C(x, x') = ( /?x /°x' ) _ ( Px Px' ) ~ ( Px Px' ) + ( PxPx' ) • This allows for a derivation 
of the kinetic-energy distribution in terms of correlation functions of vortex positions. 
Now, we can model pairing by the density-density correlation functions 

(PI W P V J A) ) = ^S(X- X ')+P^, (21) 

(pI w pT) = t^tO^-\*-A) + p^> (22) 

where Vr is the volume in which we take averages, and V\ = nX 2 is the area where the 
theta function equals one, measuring the correlation regime of vortices and antivortices. 
The contributions 



M (M — 1) 
Vr(Vr-VaY 



P^ = TTTTr rMl x - x '|- A ) (23) 



take into account that, besides pairing, vortices and antivortices keep a minimum 
distance A in the dilute gas. This is due to vortex-vortex repulsion and fast vortex- 
antivortex annihilation on small distances. The functions P x ,x' cancel out in Eq. (20). 6 
From this ansatz, two scaling regimes can be found (Nowak et ai, 2012). In the 
case of pairing, the flow field far away from the cores is given by the field of a vortex 
pair which decays as r~ 2 , and the low- momentum power law is dominated by the flow 
of random vortex pairs rik ~ k~ 2 . Above fc pa ; r ~ tt/A, the distribution exhibits the 
scaling of an ensemble of independent vortices, n^. ~ fc~ 4 , up to the healing- length 

6 If different avoidance scales A apply for vortices and antivortices, the terms do not cancel, but 
the remaining term does not alter the results for pair scaling derived here. 
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Fig. 3 Radial momentum distributions for a set of infinitely thin (£ = 0) randomly dis- 
tributed vortices. Note the log-log scale. At the highest momenta, the total number spectrum 
n(k) (dots) and the incompressible component nj(fc) (lines) follow the k~ 4 scaling of inde- 
pendent vortices. However, at momenta k < fey, with fcy(iVv = 10) ~ 1CP 2 , fcy(iVv = 100) 
~ 3 x 10~ 2 , fcy(A r v = 1000) ~ 10" 1 , n(k) becomes flat while rii(fc) continues to rise. 

scale k(: above which one can observe the vortex-core scaling <~ fc~ 6 . The above results 
show that, in a vortex dominated flow, particles with low momenta are found far away 
from the vortex cores. Particles closer to the vortex cores pick up a higher momentum. 

We remark for conciseness, as was shown by Novikov (1975), that one can obtain 
Kolmogorov 5/3 scaling from the statistics of point vortices, by choosing the density- 
density correlation functions to decay as 

^ (A) ^ (A) )-|x-xr« + Px,x', (24) 

(p V(A) p A(V) ) = pxx/ (2g) 

where the contributions P x ,x' are assumed to be equal. The integral in Eq. (20) is 
convergent for 1/2 < a < 2. This includes a = 4/3 which gives n(k) ~ fc~ 4 - 66 and 
thus Kolmogorov scaling E{k) ~ fc~ 5 / 3 . Note that the presence of vortex-antivortex 
correlations destroys the 5/3 scaling as discussed by Bradley and Anderson (2012). 

Infrared cutoff, onsider a Bose gas with density ny containing a random distribu- 
tion of vortices of either sign. We expect a decay of the coherence {(f)* (x)</>(x')) over a 

1/2 

distance |x — x'| of the order of the mean vortex distance ly = n v 1 , corresponding 
to a momentum scale fcy — tt/^v- That is because vortices appear on average at this 
distance and induce a rapid change of the phase angle if. For momenta k < fey, the 
momentum distribution n(k) needs to be sufficiently flat in order to insure convergence 
of the integral that gives the total number of particles. This restriction is not imposed 
upon the particle numbers defined by the hydrodynamic decomposition in terms of 
ni,n c ,n q . Hence, we expect a deviation of n(k) from the incompressible momentum 
distribution rii below fey. In Fig. 3, we present numerical evidence for our reasoning. 
The plot shows the single-particle (dots) as well as the incompressible momentum dis- 



14 Non-thermal fixed points: universality, topology, & turbulence in Bose gases 




Fig. 4 Illustration of Onsager's picture of thermodynamic equilibrium states of the balanced 
two-dimensional Coulomb-type vortex gas at fixed defect number. We plot the qualitative 
form of the number of available states Q(E) as a function of energy E as well as snapshots 
of typical microscopic configurations of vortices and antivortices. 

tributions (lines) for three different vortex numbers Ny. One can observe the vortex 
density dependence of the IR cutoff in the single-particle momentum distribution. The 
spectrum of the incompressible velocity field does not show this feature. Instead, the 
fc~ 4 scaling persists all the way to the lowest momenta. 

Onsager model: Non-thermal fixed point as a maximum-entropy state. We 

close with the picture Onsager (1949) developed of thermodynamic equilibrium states 
of a fixed number of vortices and antivortices in two dimensions. He used the Hamil- 
tonian of vortical flow (Lin, 1941), 

1 M 

H = -—Y t K i KM\*<-*i\)> ( 26 ) 

i>j 

to describe the dynamics of a system of M vortices in a superfluid which hence interact 
like charge carriers in a Coulomb gas. Here, the position of the i-th vortex is denoted as 
Vi = (xi, yi). Due to the fact that the x and y coordinates of each vortex are canonical 
conjugates of each other, phase space is identical with position space. Hence, for M 
vortices moving in a 2D volume V the total phase space is V . The Hamiltonian (26) 
implies that low-energy configurations feature vortices of opposite sign close to each 
other, whereas high-energy configurations require vortices of equal sign to group. Due 
to these constraints, the number of configurations Sl(E) available for the system at a 
given energy E decreases towards high and low energies, with a maximum at some 
intermediate E — Eq. This concept is illustrated in Fig. 4. According to Boltzmann, 
the entropy is 

S(E)=\n[n(E)} , (27) 
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and the inverse temperature 1/T = dS/dE is positive for E < E and negative for E > 
E . It follows that positive-temperature states are characterised by vortex-antivortex 
pairing, while negative-temperature states feature vortices of the same circulation to 
cluster. At the point of maximum entropy S(Eq) and infinite temperature, Onsager 
expected a state of uncorrelated vortices and antivortices. 

3-dimensional systems. We briefly comment on the case of vortex lines and loops 
in three dimensions. A formulation similar to the Onsager point vortex model is pos- 
sible (Nemirovskii, 1998; Nemirovskii et al, 2002; Tsubota, 2008). 

For taking into account the most general case of squeezed vortex loops, it is helpful 
to consider elliptical filaments, characterised by a major radius r a and minor radius rf>. 
Three scaling regimes can be distinguished. For the lowest momenta, one has n(k) <~ 
fc~ 2 , which equals the infrared scaling in the presence of a vortex ring. For momenta 
k a <C k -C fc&, one finds n(k) <~ fc -3 , which coincides with the infrared scaling of 
two anti-circulating vortex lines. For the ellipse, above kb, the momentum distribution 
scales like n(k) <~ fc~ 5 . This is the scaling of a single vortex line and can also be found 
as the high-momentum scaling of a vortex ring or a pair of straight vortex lines. For 
more details we refer to (Nowak et al, 2012). 



3 Non-thermal fixed point of a vortex gas 

In extending the concept of universality to time evolution far from thermal equilibrium 
one expects that also away from the thermal limit the character of dynamical evolu- 
tion can become independent of the microscopic details. Looking at closed systems this 
implies that, in approaching critical configurations, the evolution must become inde- 
pendent of the particular initial state the system has started from and critical slowing 
down in the actual time evolution is observed. Time evolution near the fixed point 
becomes equivalent to a coarsening transformation. It should look like a self-similar 
pattern which is observed through a microscope while one continuously turns the mag- 
nification of the lens to smaller focal lengths. Precisely at the fixed point the system 
becomes stationary due to its self-similarity under time translations. Considering a 
generic isolated system, it can evolve into the vicinity of a non-thermal fixed point, 
stay there over long times before it eventually undergoes thermalisation. This picture 
is supported by the derivation of the strong turbulence scaling laws in the frame of 
renormalisation-group theory for correlation functions (Berges and Hoffmcister, 2009). 

In this section we discuss mostly numerical results obtained for the time evolution 
of an isolated two-dimensional gas. This is characterized, from the turbulence point 
of view, by direct and inverse cascades, with fluxes determined by local conservation 
laws in Fourier space. From the perspective of defect formation a diluting ensemble 
of vortices and antivortices marks the approach of the non-thermal fixed point. We 
identify the mechanism for this dilution process and demonstrate that it is consistent 
with the coarsening transformation picture of universal dynamics. By reducing the 
characterization of the momcntaneous configuration to a few length parameters we can 
illustrate the slowing of the evolution near the fixed point. We show that it depends 
on the chosen initial state how closely the critical point is approached. 
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Fig. 5 Phase angle (p(x, t) (upper panel) and spatial density n(x) (lower panel) at six times 
during a single run of the simulations in d — 2 on a space-time lattice with side lengths 
L = N s a 3 , lattice spacing a 3 , imposing periodic boundary conditions. The GPE is written in 
terms of the dimensionless variables: g = 2mga 2 a ~ d , t = t/{2ma1) and <p(t) = ^ya|exp(2it). 
Parameters are: q = 3 ■ 1CT 5 , N = 10 s , N s = 256. Shown times are (by line, each from left to 
right) : t = 0: Ordered phase at initial preparation, t = 104: During build-up of strong density 
and phase gradients, t = 210: Formation of vortices and antivortices. t = 420: Unbinding of 
vortex-antivortex pairs, t — 6550: Slowing down of vortex dynamics near the non-thermal 
fixed point, t — 10 s : Few seperated vortex antivortex pairs close to non-thermal fixed point. 



3.1 Time evolution of vortex patterns and momentum spectra 



In this section we have a closer look at the process of vortex formation and of the 
Bose gas approaching the non-thermal fixed point, concentrating again on the two- 
dimensional case. Vortical excitations can be created in large numbers, e.g., within 
shock waves forming during the non-linear evolution of a coherent matter-wave field. 
We follow the exemplary evolution of phase and density profiles in Fig. 5. Six snapshots 
are shown, taken at the dimensionless times i as indicated in the caption. The initial 
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Fig. 6 Single-particle, angular and ensemble averaged mode occupation numbers 
n (k) = / d d ~ 1 Qk (4>* (k)0(k)) cnsom bic as functions of the radial momentum k, at 
the same evolution times as in Fig. 5. The dimensionless lattice momenta are 
k = EiLi 4sin 2 (fei/2)] 1//2 , ki = 2niii/N 3 , iii = —N 3 /2,...,N 3 /2. Parameters are the same 
as in Fig. 5. Note the double-logarithmic scale. At early times (top left), scattering between 
macroscopically occupied modes leads to excitation transfer to high momenta. Vortex forma- 
tion (t ~ 210) sets in when this process has reached the healing length scale fej ~ 1.26. At 
t > 10 4 a bimodal power law emerges, characteristic for the non-thermal fixed point. 

field configurations were prepared by macroscopically populating a few of the lowest 
momentum modes in the computation such that the resulting condensate density in 
position space varied between zero and some maximum value. One observes strong 
phase gradients forming due to the non-linear evolution. At around t ~ 100, these gra- 
dients produce shock fronts delimited by phase defects which in the following collapse 
into vortex trains. Scattering processes between vortices quickly isotropise phase and 
density fluctuations. 

In Fig. 6, we show the corresponding time evolution of the angle- and ensemble- 
averaged radial momentum spectrum. The early times (top left) are characterised by 
scattering between macroscopically occupied modes. Once excitations of the order of 
the healing length £ are created, vortex formation sets in. Shortly after vortices are 
created the spectrum exhibits a power-law behaviour 2.85 < £ < 3.0 within a range of 
momenta k G [0.04 : 0.4], see the lower left panel of Fig. 6. Subsequently, the evolution 
slows down and a quasi-stationary period is entered. During an intermediate stage 
(bottom centre and right panels of Fig. 6) of the vortex-bearing phase two distinct 
power laws develop which are in excellent agreement with the analytical prediction in 
Eqs. (5) and (7). While in the ultraviolet the exponent Cp V = d exhibits weak wave 
turbulence, Eq. (5), in the infrared, the exponent confirms the field theory prediction 



IB 

Q 



d + 2, cf. Eq. (7). During the ensuing evolution, the weak-wave-turbulence 



scaling decays towards £ = 2, reflecting a thermal UV tail. Note that in d = 2, the 
weak-turbulence exponent £p V = 2 is identical to that in thermal equilibrium in the 
Rayleigh- Jeans regime, n(k) ~ T/k 2 (Zakharov et ai, 1992). In d = 3 we observe, at 
late times, a change of the infrared scaling behaviour from £ = <i + 2 = 5to£ = 3, 
pointing to the development of pairing correlations (Nowak et at, 2012). 
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Fig. 7 Direct kinetic-energy and inverse particle fluxes in d — 2 at t — 6550 (Nowak el al., 
2012). Parameters are as in Fig. 6. Note the logarithmic fc-axis. A positive kinetic energy flux 
is seen in the UV, a negative particle flux in the IR. Units: [P] = [gnQ] — (4m 2 af +4 ) _1 . 

At late times, after the last vortical excitations have disappeared, we observe the 
entire spectrum to become thermal, i.e., exhibit Rayleigh- Jeans scaling with £ = 2 (not 
shown). We emphasise that thermal scaling of the single-particle occupation number 
has C — 2 despite the fact that quasi-particles with a linear dispersion are expected to 
thermalise in the regime of wave numbers smaller than the inverse healing length. 

3.2 Local transport in momentum space and inverse particle cascade 

From the above findings the question arises, why the system selects the particular 
exponents £p V = d and (q 1 = d + 2 from the set of four possible exponents given in 
Eqs. (5) and (7). For this, the fluxes underlying the stationary but non-equilibrium 
distributions are relevant (Berges and Sexty, 2011). The timeline of distributions shown 
in Fig. 6 suggests that the evolution of the gas involves transport of particles originating 
from the intermediate momentum regime k ~ 0.05 . . . 0.2, which during the initial 
evolution gets strongly overpopulated. The particles drift both towards lower and higher 
wave numbers, building up a bimodal power-law distribution. To describe the character 
of this bidirectional flux we plot, in Fig. 7, the radial particle and kinetic-energy flux 
distributions Qk and Pk , respectively, at t — 6550 corresponding to the bottom center 
panel of Fig. 6. Note that the radial particle flux density Qk is multiplied by gn to 
have the same units as the energy flux density P(k). These flux densities are defined 
through the balance equations (3) and (4), respectively, with kinetic energy density 
£fc = rikk 2 /2m. They are determined by integrating the numerically obtained particle 
and energy spectra N(k) and E(k) up to the scale k. 

The graph supports the interpretation of the transport in terms of an inverse 
particle cascade in the IR and a direct energy cascade in the UV (Nowak et at, 2012), 
in accordance with the appearance of the bimodal momentum distributions in Fig. 6. 
Although the derivation of the IR exponents requires the full dynamical theory with 
non-perturbatively resummed self-energies, the signs of the fluxes correspond to the 
respective scaling exponents, i.e., (q 1 in the IR, and £p V in the UV. Moreover, at late 
times, the kinetic-energy flux P almost vanishes due to a thermalised UV momentum 
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Fig. 8 (a) Vortex density p as a function of time Evolution for various initial conditions 
(Schole et al, 2012), averaged over 20 runs on a grid of size N s — 1024. Lines show different 
power-law evolutions. The vortex density follows power laws p(t) ~ t~ ai with two different 
exponents cti, i — 1,2. The exponent during the early stage depends considerably on initial 
conditions 1 < a\ < 2, whereas the late stage features a decay exponent in a narrow interval 
0.3 < «2 % 0.4. The closest approach to the non-thermal fixed point is at t ~ (5 . . . 10) x 10 . 
(b) Sketch of different scattering events between two vortex-antivortex pairs in d — 2 di- 
mensions. Left panel: Scattering of two vortex pairs resulting in mutual annihilation of two 
vortices and the emission of density waves. Right panel: Scattering of two vortex pairs, lead- 
ing to a change in the vortex-antivortex distance Id and pair velocity. We remark, that once 
Id ~ £ a vortex pair decays rapidly under the emission of density waves. 

distribution, but Q still reshuffles particles and therefore energy, with the zero mode 
acting as a sink, keeping the system out of equilibrium close to the non-thermal fixed 
point. 

We finally remark that a necessary condition for a non-equilibrium stationary dis- 
tribution is energy damping at large k (Zakharov et al., 1992). Moreover, energy and 
particle number conservation in the interaction of different momentum modes can be 
shown to imply the existence of at least one more sink, i.e., a region where the right- 
hand sides of Eq. (3) effectively has an additional damping term ~ T(k)n(k), with 
negative T. In between these sinks, a source region supplies the input to the bidirec- 
tional flux pattern towards the UV and IR. Using kinetic theory, one can show that 
under certain conditions a positive /c-independent flux transports energy, P > 0, while 
a negative flux transfers particles, Q < (Zakharov et at, 1992), Remarkably, this 
pattern remains valid in our case, besides the UV weak- wave-turbulence regime also in 
the IR region where the exponent emerges from a fixed point of the non-perturbative 
dynamic equations for Green's functions. As already pointed out by Scheppach et al. 
(2010), the derivation of the IR exponent = d + 2 requires sufficiently well defined 
quasi particles, suggesting a treatment in terms of the Quantum Boltzmann equation 
with a momentum dependent scattering matrix element to be applicable. From this 
point of view, the negative flux Q and scaling in the IR and the positive flux P and 
weak wave turbulence in the UV, as observed in the numerics, emerge as a necessary 
consequence of conservation laws and transport processes described by wave-kinetic 
transport equations with non-trivial interactions. 
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Fig. 9 Normalised pairing correlation functions gvA, see Eq. (28), as a function of the vor- 
tex-antivortex distance r, for six different times ti (in lattice units), (a) <?va(t") during the 
rapid-decay stage, averaged over 174 runs on a grid of N s — 1024. Inset: Vortex density p 
time evolution, from Fig. 8a. (b) <;va(»~) during the slow-decay stage, same averaging as in 
(a). Inset: Evolution of mean vortex-antivortex pair distance In- 

3.3 Approach of the non-thermal fixed point and critical slowing down 

Let us study in some more detail the evolution of the system towards and away from 
the non-thermal fixed point and focus on universal aspects of the dynamics of vortices. 
For details we refer to (Schole et ai, 2012). Fig. 8a shows the time evolution of the 
vortex density p(t) = (N Y (t) + N A (t))/V, where N y( - A 1(t) is the mean number of 
vortices (antivortices) in the volume V at time t, for different specific choices of the 
initial state, cf. (Schole et ai, 2012). In all runs, vortex formation occurs around 
ty ~ 10 3 , apparent from the steep increase of vortex density around this time. For 
t > ty, two distinct stages in the vortex density decay are observed, a rapid early 
stage and a slow late stage. We have repeated our simulations on various grid sizes, 
N s e [256, 4096]. Thereby, we found that decay exponents saturate for and above 
N s = 512. We attribute deviations on smaller grids to effects from regular (integrable) 
dynamics of few- vortex systems (Aref, 1983). We remark that the onset of the slow 
decay coincides with the development of a particular scaling behaviour in the single- 
particle momentum distribution n(k) ~ fc -4 , which by Nowak et al. (2011, 2012) was 
shown to signal the approach of the non-thermal fixed point and the formation of a set 
of randomly distributed vortices. In this context, the reduction of the vortex density 
decay exponent, compared to the early stage of rapid decay, is interpreted as due to 
(critical) slowing down of the nonlinear dynamics near the non-thermal fixed point. 

We can discuss the dynamical transition in the vortex annihilation dynamics in 
terms of characteristic features of the vortex-antivortex correlation function 



5va( 



(p v (x J t)p A (x',t)) 



(28) 



where p v ( A '(x,£) = E^( x 



(pV (x , t))((0 A (x , )t)) ' 

x.j(i)) is the distribution of (anti) vortices in a single 
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Fig. 10 Trajectories of multi- vortex states in the space of inverse coherence length 1/lc 
and inverse mean vortex-antivortex distance 1/Zd, starting from t — ty. Comparison of the 
thermal line ^(Ty^fT)) for a range of temperatures T (right) with the corresponding 
dynamical evolution (left). Dashed lines mark the minimal values 2/L = 0.018£ , available 
on a grid of size N s = 512. Note that the (l/foj-axis interval [0.25,0.55] has been cut out. 

run at time t. For sufficiently large ensembles, <7va is a function of r = |x — x'| 
only. At early times, one finds a strong pairing peak in gvA^jt) near r = 0, see 
Fig. 9a. This peak gets quickly reduced and a hole is 'burned' into the correlation 
function near the origin, see Fig. 9b. Following the time evolution of the spatial vortex 
distribution we observe that this involves qualitatively different processes: Mutual 
annihilations of closely positioned vortices and antivortices occur under the emission 
of sound waves. Further separated vortices can approach each other in different ways as 
illustrated in Fig. 8b. The scattering of two pairs can directly lead to the annihilation 
of one pair under the emission of sound waves. This includes events where the dipole 
length reduces below a certain threshold, implying a density dip rather than a vortex- 
pair. This density dip can still interact with other vortices but will quickly vanish. 
Alternatively, the scattering reduces the vortex-antivortex separation within one pair 
while it increases it within the other, in accordance with the Onsager point-vortex 
model (Onsager, 1949). We refer to this characteristic change in flVA^) as a vortex 
unbinding process. Recall that the Onsager model does not contain a kinetic term for 
the point defects such that their interaction potential energy can not induce relative 
acceleration of the vortices. Changes in their separation must occur dynamically via 
collisions. 

At around the time t% < t < t± the power-law exponent of the vortex density decay 
changes to about a third of its previous value, see the inset of Fig. 9a. Computing 
the mean vortex-antivortex pair distance /rj , by averaging over distances between each 
vortex and its nearest antivortex, we find that, in accordance to the previous discussion, 
lj) grows continuously, exhibiting two characteristic stages, see the inset of Fig. 9b. At 
times t > 10 4 , lr>(t) approaches the power-law solution Z D ~ p -1 / 2 , as expected for 
uncorrelated vortices. 
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We close by focusing on the growth of long-range coherence, associated with the 
annihilation of topological defects (Levich and Yakhot, 1978; Kagan et ai, 1992; Ka- 
gan and Svistunov, 1994; Damle et ai, 1996; Berloff and Svistunov, 2002; Svistunov, 
2001; Kozik and Svistunov, 2009; Nazarenko and Onorato, 2006). From this point of 
view, freely decaying superfluid turbulence is a particular example of phase-ordering 
dynamics after a quench into the ordered phase (Bray, 1994). Whereas in three dimen- 
sions, a second-order phase transition connects a normal-fluid and a superfluid phase, 
a Bose gas in two dimensions experiences a Berezinskii-Kosterlitz-Thouless (BKT) 
transition (Berezinskii, 1971; Kosterlitz and Thouless, 1973). For the two-dimensional 
ultracold Bose gas, experimental and theoretical results support the understanding of 
the phase transition in terms of vortices undergoing an unbinding-binding transition 
(Hadzibabic et ai, 2006; Simula and Blakie, 2006; Schweikhard et ai, 2007; Giorgetti 
et al, 2007; Weilcr et al, 2008; Bisset et al, 2009; Foster et ai, 2010). 

In this context, we are interested in a comparison between correlation properties 
observed in the non-equilibrium dynamics near a non-thermal fixed point and those 
known from equilibrium studies. We compute the dynamical trajectory of the vortex 
gas in the space of inverse coherence length and inverse mean vortex-antivortex pair 
distance. We compare our results to simulations of a thermal two-dimensional Bose 
gas specifically for our system parameters. We define a coherence length Z£ in terms 
of the integral over the angle-averaged first-order coherence function, 
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Iq measures the spatial extension of the first-order coherence function. Contrary to 
?"coh = / dr r 2 g^\r) / J drrg^\r), the quantity Iq does not sum up values of g^\r) 
weighted by the distance, which would enlarge insignificant contributions at large r. 
In addition, it does not overestimate the coherence of flat distributions. In equilibrium 
this quantity smoothly interpolates between the regime of exponential decay of ~ 
exp(— r/£c) above the BKT transition and its power-law decay below. 

In Fig. 10, we follow the time evolution of the gas for t > ty. One can observe 
that a state of low coherence and small mean vortex-antivortex pair distance evolves 
towards larger coherence and larger vortex-antivortex separation. As discussed above, 
this is due to vortex annihilations and vortex-antivortex unbinding. For times t > 10 4 , 
the coherence length grows as Iq ~ p -1 / 2 , in the same way as Id shown in Fig. 9b. The 
evolution considerably slows down for ~ 1/Zd — >• 0, when the gas starts to show 
the characteristic scaling n(k) ~ fc -4 . After spending a long time near the non-thermal 
fixed point, Id declines because the last remaining vortex-antivortex pairs reduce their 
size prior to their annihilation and the equilibration of the system. Our understanding 
of the non-thermal fixed point as a configuration with a few, maximally separated 
pairs on an otherwise maximally coherent background implies it to be located near 
the crossing of the dashed lines. Hence, the non-thermal fixed point is approached most 
closely between t~5x 10 5 and t ~ 10 6 — t 3 . To set the above evolution in relation 
to equilibrium configurations, we also show the thermal line {Zj^ (T), (T)} for a 
range of temperatures T for which the zero-mode population does not vanish. 

Our results allow to draw a picture of the evolution path towards and away from 
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Fig. 11 Momentum distributions of a 3-dimensional gas for different initial conditions. Left: 
Initial conditions with varying high momentum decay k~ a , characterised by the paramter 
a — 2.5, 3, 4, 6, 10, representing different evaporative cooling quenches. Right: At intermediate 
times, when the system is closest to the non-thermal fixed point, showing the steep 1R power 
laws ~ k for initial conditions representing a strong cooling quench, a > 3. 

the fixed point. The non-thermal fixed point is characterised by a few pairs - in the 
extreme case one pair - of far-separated anti-circulating vortices and bears similarities 
with the equilibrium BKT fixed point. However, while the phase transition also features 
unbinding of vortices, finite temperature implies the simultaneous excitation of many 
rotons, i.e. strongly bound vortex-antivortex pairs. The non-thermal fixed point is 
clearly identified by strong wave turbulent scaling in the infrared limit, n{k) ~ fc -4 . 
The high-energy modes are much weaker populated, e.g., far below the BKT critical 
temperature or remain out of equilibrium. The details of the UV mode populations 
are determined by the way the non-thermal fixed point is being approached. 

The way the system is forced, here, to approach the non-thermal fixed point gener- 
alises the protocol of Kibble (1976) and Zurck (1985). A strong sudden quench replaces 
the adiabatic approach of the BKT transition. The dynamical evolution in the vicinity 
of the BKT critical point was studied by Mathey and Polkovnikov (2009), in terms of a 
perturbative renormalisation group analysis. The route to a non-perturbative analysis 
in the strong-coupling regime is provided by out-of-equilibrium functional techniques, 
see Refs. given at the end of Sect. 2.2. 

3.4 Dependence on driving and a new route to Bose condensation 

Let us finally give a taste of the relevance of the strength with which the system is 
driven away from thermal equilibrium for the approach of the non-thermal fixed point 
and, as a by-product, find new aspects of how superfluid turbulence affects the process 
of Bose condensation (Nowak and Gasenzer, 2012). In this section we restrict ourselves 
to a gas in d — 3 dimensions. We choose the overpopulated momentum distribution 
of the gas directly, eliminating the early instability phase of our previous simulations 
during which overpopulation is induced through non-linear scattering of wave modes. 
The initial field in momentum space, </>(k, 0) = \J n(k, 0)exp{i<y9(k, 0)}, is parametrized 
in terms of a randomly chosen phase ip(k, 0) G [0, 2ir) and a density n(k, 0) = /(fc)i^, 
with i/k > drawn from an exponential distribution P(vk) = exp(— v-^) for each k. The 
spectrum is flat at low k and falls off according to the function f(k) = / a /(&o + ^ a )> 
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Fig. 12 Depending on the strength a of the initial cooling quench, the gas can thermalise in a 
near-adiabatic manner to a Bose-Einstein condensate. Alternatively, it can first approach and 
be critically slowed down near a non-thermal fixed point where it is characterised by a scaling 
spectrum n(k) ~ fc -5 in the IR. Furthermore, dynamical scale separation of the incompressible 
(red points) and incompressible (blue points) components of the gas (cf. Nowak et al. (2012), 
and last par. of Sect. 2.3) occurs. The incompressible component corresponds to transverse 
vortical flow. The grey points indicate the quantum-pressure component. 

where a controls the deviation from a thermal decay with a — 2. We choose a cutoff 
(koO a — 0.2/0.44" and normalization f a — 400/0.44". We compare results for a range 
of different cooling quenches denned by the power-laws a = 2.5, . . . , 1 0.0, varying the 
total number between N = 10 9 (a = 2.5) and N = 4.3 x 10 8 (a = 10). 

In Fig. 11 we show n(k,t) over the radial momenta k = |k| at the initial time as 
well as at a moment when the system is closest to the non-thermal fixed point, for the 
different choices of a. During the initial evolution (t < 10 2 r) the mode occupations 
gradually spread to lower wave numbers, at the same time depositing energy into the 
high-momentum tail. We emphasise that cutting away sufficiently much population 
at high momenta initially is necessary if the system is supposed to approach the 
non-thermal fixed point during its rethermalisation: As described in Sect. 3.2, the 
approach of the fixed point is characterized by a dual cascade in which the energy of 
the intermediate-fc overpopulation gets deposited in the high-fc tail, carried there by 
a few particles, while the majority of overpopulation particles moves towards the IR, 
conserving overall energy and particle number. Hence, only strong cooling quenches 
allow for the build-up of a steep population far into the IR. 

At late times, the spectra developing from the different initial a differ strongly. 
For a > 3, the distribution develops a bimodal structure, with a power-law behaviour 
n(k) ~ fc~ 5 in the infrared (IR) and n(k) ~ fc -2 in the UV. At very long times, 
this bimodal structure decays towards a global n(k) ~ k~ 2 (not shown). For a < 
3, the distribution directly reaches a thermal Rayleigh- Jeans scaling n(k) ~ T/k 2 . 
Preliminary results (Nowak et at, 2013) show that the trajectories, when plotted as 
in Fig. 10, approach the non-thermal fixed point in the lower left corner the closer, 
the larger a is chosen. While a = 2.33 leads to a trajectory near the thermal states, 
a = 10 induces a motion similar to the one of the black points. 
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Fig. 13 Momentum distributions for N s randomly distributed black, i.e., static solitons 
(symbols), showing a k~ 2 power law behaviour between the scale k„ s given by the mean 
soliton separation and kkink given by the healing length. The spectra are found to be in 
excellent agreement with analytically predicted lines (Schmidt et al, 2012). 

Much experimental effort is undertaken at present to study the dynamics of conden- 
sation under conditions of rapid evaporative cooling (Pitaevskii and Stringari, 2003; 
Ritter et al, 2007; Weiler et al, 2008; Smith et at, 2012). Moreover, many-body 
dynamics of coherent bosonic excitations is intensively studied in solid-state systems 
consisting of magnons (Demidov et al, 2007; Demidov et al, 2008; Keeling and Bcrloff, 
2008; Nowik-Boltyk et al, 2012) or polaritons (Saba et al, 2001; Kasprzak et al, 2006; 
Balili et al, 2007; Lagoudakis et al, 2008; Amo et al, 2011; Keeling and Berloff, 2011). 
Recently, condensation has been discussed for the case of gluons as an intermediate 
stage of heavy-ion collisions (Blaizot et al, 2012; Berges et al, 2012) and, for relativis- 
tic scalars, been found to rely on a nonperturbative inverse particle cascade (Berges 
and Sexty, 2012). 

Our results confirm that Bosc-Einstein condensation in a non-equilibrated and 
under-cooled gas can have the characteristics of a turbulent inverse cascade (Svis- 
tunov, 1991; Kagan et al, 1992; Kagan and Svistunov, 1994; Berloff and Svistunov, 
2002), corresponding to a quasi- local transport process in momentum space, into the 
low-energy modes of the Bose field. The main new finding is that the supcrfluid turbu- 
lence period can appear in two different forms (Nowak and Gasenzer, 2012). The two 
possible paths to Bose-Einstein condensation are shown schematically in Fig. 12. If a 
sufficiently small amount of energy is removed, a thermal Rayleigh- Jeans distribution 
forms in a quasi-adiabatic way. The chemical potential increases, and a fraction of 
particles is deposited in the lowest mode, forming a Bose-Einstein condensate. In the 
second scenario, after a sufficiently strong cooling quench, the system develops tran- 
sient scaling behaviour in the momentum distribution prior to condensate formation, 
i.e., it approaches the non-thermal fixed point. These scenarii differ qualitatively in 
how the condensate mode builds up as a function of time (Nowak and Gasenzer, 2012). 
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Fig. 14 Time intervals of a single run of the classical field equation (10), showing solitons 
which oscillate inside a trapped one-dimensional ultracold Bose gas. The gas is initially non- 
interacting and thermalised, with T — 360a;h o , in a trap with oscillator length = 8.5 (in 
grid units). At time t = the interaction is switched to gm = 7.3 x 10~ 3 , and a cooling 
period using a high-energy knife is applied. The panels show the one-dimensional colour-en- 
coded density distribution as a function of time. Left panel: The initially imposed interaction 
quench causes strong breathing-like oscillations and the creation of many solitons. Right panel: 
Breathing oscillations have damped out, leaving a dipolar oscillation of the bulk distribution, 
with clearly distinct solitons in the trap. For details cf. Schmidt et al. (2012). 

4 Other systems 

4.1 Soliton ensembles 

In the preceding sections, it was discussed how strong wave turbulent scaling of the 
momentum spectrum of a degenerate Bose gas can be understood from the statistics 
of vortices. In this context, correlations between vortices and antivortices play a cru- 
cial role. In three spatial dimensions similar relations exist between scaling and the 
creation of vortex lines and rings (Nowak et al., 2012). Looking at systems in one 
spatial dimension we find that solitons, in particular dark solitons play a crucial role 
in realising strong IR wave turbulence. Solitons interact with other defects as well 
as sound excitations and show characteristics of localised quasi-particle excitations. In 
the following we briefly summarise to what extent random ensembles of solitons can be 
seen as characterising a non-thermal fixed point in one-dimensional bosonic systems. 

A model of randomly positioned grey solitons, being solutions of the Gross-Pitaev- 
skii equation (10) in the vicinity of each soliton density dip, was discussed by Schmidt 
et al. (2012). In the limit of large soliton separations compared to the healing length, 
analytic expressions for the momentum spectra were derived in a homogeneous sys- 
tem as well as under the constraint of a trapping potential. A central result is de- 
picted in Fig. 13 for the case of randomly distributed black solitons in a homogeneous 
background. The momentum distributions show a k~ 2 power law behaviour between 
the scale k ns corresponding to the mean soliton separation and fckink which marks 
the healing length scale. At higher momenta the distribution features an exponential 
decay following from the characteristic way the spatial density drops inside the soli- 
ton cores while at low momenta a flat distribution marks the long-range exponential 
decay of the coherence. The soliton configurations in one dimension resemble the tur- 
bulent phenomena in d = 2 and 3 dimensions. Their presence is accompanied by a 
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quasi-stationary power-law momentum distribution which marks the self-similarity of 
a random distribution of sharp phase kinks: looking at the system within a window, 
the size of which is below the mean soliton separation and well above the core width, 
the system looks the same irrespective of the window size. There is either a kink seen 
within the window or not. 

Transport in momentum space, in the form of an inverse cascade towards the IR, 
corresponds to soliton decay which increases the inter-soliton distance and extends 
the self-similar regime towards smaller momenta (see Fig. 13). We remark that the 
predicted strong wave turbulent scaling ~ k~ d ~ 2 which was found consistent with 
the vortex scaling in d = 2 and 3 does not give the scaling ~ k~ 2 for the solitons in 
d = 1. The reason for this is expected to be similar as for the case of domain walls in 
d — 2, 3, in a multi-component system as described in Sect. 4.2 below: new transport 
equations must be set up in which a different conserved current, relating to transport 
of spin wave excitations, leads to a different IR scaling. A soliton in d = 1 does not 
exhibit transverse (incompressible) superfluid flow as that around a vortex core which 
decays as |v(r)| <~ 1/r. In this context we remark that the longitudinal (compressible) 
component in d = 2 and 3 dimensions (Fig. 12) shows a by one weaker IR exponent. 

The formation and evolution of soliton excitations in trapped one-dimensional Bosc 
gases can be studied by means of the classical field equation. A possible scenario of 
far-from-equilibrium dynamics involves an initially non-interacting thermal gas that is 
quenched by a sudden ramp of the interaction, as studied by Schmidt et al. (2012). To 
allow the emerging collective excitations to form solitons at a desired density, evap- 
orative cooling helps to achieve the required reduction of the UV mode populations. 
An exemplary position-space evolution is depicted in Fig. 14. 

4.2 Domains and defects in two-component systems 

As a further example we briefly comment on what is already known in the con- 
text of strong wave turbulence and defect formation in multi-component, e.g., spinor 
gases (Karl et al. (2013); Fig. 15). A two-component Bose gas with contact interac- 
tions is well known to possess two different ground states depending on the value of 
the parameter a — 511322/912 which parametrizes the relative strength of intra-species 
couplings gu and 922 as compared to inter-species interaction g\ 2 (Timmcrmans, 1998; 
Kasamatsu and Tsubota, 2006) and is the subject of increasing experimental inves- 
tigations (Sadler et al., 2006; Vengalattore et al., 2008; Guzman et al., 2011; Nicklas 
et at, 2011). For a > 1, which is called the immiscible regime, the inter-species in- 
teraction energy is greater than that of the interactions within one species. Hence, in 
the ground state of the system, the spatial overlap of the components is minimized 
through domain formation and spatial separation of the two components. In contrary, 
for a < 1, the two components become miscible and uniformly distributed when in 
the ground state. One can use the parameter a to change the properties of the system 
in the yet unexplored region of non-equilibrium quasi-stationary states. 

Dynamical instabilities serve to drive the system far from equilibrium. This leads 
to momentum distributions of the different components which are characterized by 
a strong overpopulation at intermediate momenta as compared to thermal equilib- 
rium distributions with the same energy and particle numbers. Similar to our finding 
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Fig. 15 Non-thermal fixed point and domain formation in an immiscible two-component 
Bose gas. Left: Occupation number spectrum in d = 2 spatial dimensions at different times 
(note double-log scale). The spectrum at late times shows a bimodal structure with a thermal 
scaling n(k) ~ k~ 2 at high momenta and an IR power-law n(k) ~ fc~ 3 ' 5 , similar to scaling seen 
when charge separation is found in the relativistic simulations, see Fig. 16. Centre: Polariza- 
tion {S z (x)) , i.e., density difference of the two gas components, S z (x) = (ni — ri2)/(rii + «2). 
ft demonstrates the formation of spin domains as a result of dynamic demixing. Right: Spatial 
distribution of incompressible kinetic energy, £j(x), showing high amounts of transverse hy- 
drodynamic, i.e., vortical flow around spin domains and especially around point- like defects. 

for the case of a one-component Bose gas in d = 1,2,3 dimensions, such a far-from- 
equilibrium distribution of particle momenta induces a redistribution of particles in 
momentum space both towards lower and higher momenta. Subsequently, the system 
encounters long-lived transient states with non-topological and quasi-topological de- 
fects including domain walls, vortices in a single species, and skyrmions in the coupled 
spin system. Distinguishable types of defects are produced for different values of the 
external parameter such that this allows to induce a transition between meta-stable 
non-equilibrium ordered states. In summary, one obtains an example of how to extend 
the concept of a phase transition into the realm of far-from-equilibrium time-evolution. 

4.3 Charge separation in reheating after cosmological inflation 

In cosmological models of the universe, reheating describes the epoch starting at the 
end of inflation (Allahverdi et al., 2010). During this period the potential energy of the 
inflaton field is redistributed into a homogeneous and isotropic hot plasma of particle 
excitations. These become a substantial part of the further expanding universe. Simple 
models describing reheating after inflation invoke self-interacting scalar fields. One of 
the popular scenarios involves the parametrically resonant amplification of quantum 
fluctuations of the macroscopically oscillating inflaton field. The amplified modes rep- 
resent the emerging matter content of the universe (Kofman et al., 1994; Traschen and 
Brandenberger, 1990). Various theoretical approaches have been proposed to model re- 
heating. As both, the inflaton and the amplified modes are strongly populated, classical 
field simulations can be applied to describe their evolution (Khlebnikov and Tkachev, 
1996; Prokopec and Roos, 1997; Tkachev et al, 1998). 

We focus on a scenario of parametric resonance in a globally 0(2) or, equivalently, 
£7(1) symmetric relativistic scalar field theory (mass parameter m = 0) in d = 2, 3 di- 
mensions (Gasenzer et al, 2012). Shortly after the resonant excitations have set in, a 
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Fig. 16 Non-thermal fixed point and charge separation in the nonlinear Klein-Gordon equa- 
tion for a complex field. Left: Occupation number spectrum in d — 2 spatial dimensions at 
different times. The double-log scale exhibits the bimodal power-law, with n(k) ~ fc -3 ' 5 at 
low momenta. Center: Modulus of the field, </>(x) showing worm-like regions (dark grey) of 
near zero field. Right: Charge distribution jo(x) of the field, showing homogeneous regions of 
opposite charge. 

spatial separation of charges occurs as shown in Fig. 16 (right panel). Both, charge and 
anti-charge overdensities become uniformly distributed within slowly varying regions 
which are separated by sharp boundary walls of grossly invariant thickness. These 
walls have a character similar to topological defects and appear for initial conditions 
corresponding to the parametric reheating scenario. One observes that the presence 
of charge domains coincides with the visibility of non-thermal stationary scaling so- 
lutions, see Fig. 16 (left panel), as discussed before within the context of ultracold 
atomic gases. In this way, a link is established between wave turbulence phenomena 
as discussed, e.g., by Micha and Tkachev 2003; Micha and Tkachev 2004; Berges et al. 
2008; Berges and Hoffmeister 2009; Scheppach et al. 2010; Berges and Sexty 2011; 
Carrington and Rebhan 2011, and long-lived quasi-topological structures in the infla- 
ton field. Strong non-thermal stationary scaling solutions have also been observed for 
the case of 0(4)- and O(10)-symmetric scalar fields in d = 3,4 dimensions (Berges 
et al, 2008; Berges and Sexty, 2011). The latter work raises the interesting question 
after the corresponding spatial configuration at the fixed point. 

5 Outlook 

In these notes we have pointed out the possibility of a universal duality between decay- 
ing topological defects and a non-perturbative inverse wave-turbulent cascade. This 
cascade requires the generation of (quasi-)topological configurations far from thermal 
equilibrium and their slow decay, going together with an increase of coherence and 
defect separation. Under these conditions, we expect power-law scaling in a regime 
between the scales l/£, £ being the microscopic extent of the defect core, and 1/Zcoh, 
where Icon is the coherence length approximated by the mean distance between de- 
fects. In this setting, an inverse particle cascade is generated by defect dilution, and the 
associated power-laws can be found from the scaling properties of the respective single 
defect. We have shown this mechanism to exist in soliton- and vortex dominated single- 
component Bose gases (Nowak et al., 2011; Nowak et at, 2012; Schmidt et al., 2012; 
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Table 1 Strong wave turbulence scaling of momentum spectra of different distributions as 
expected at a non-thermal fixed point, for different types of defects (first column) in various 
systems (second column) in d dimensions. Concerning the definitions of the defects we refer 
to (Rajaraman, 1982; Kasamatsu et al, 2005). The 'Skrymion' denotes the respective defect 
arising in the non-linear sigma model (Anderson- Toulouse vortex). 

Nowak and Gasenzer, 2012), decaying domain walls and vortices in two-component 
Bose gases (Karl et al., 2013), as well as transient charge domains in complex rela- 
tivistic scalar theory (Gasenzer et al., 2012). It is emphasised that the stability of these 
defects does not need to be topological, as the examples of solitons and charge domains 
indicate. Cf., e.g., Lee and Pang (1992) for a review of non-topological solitons. A list 
of expected scaling properties is given in Table 1. 

A variety of (quasi-)topological excitations are known to exist in superconductors, 
magnets, and cosmic fields (Lee and Pang, 1992; Thouless, 1998; Pismen, 1999; Nel- 
son, 2002). Specific examples are monopoles in gauge fields (Rajantie, 2002; 2012) 
and exotic magnets (Castelnovo et al., 2008), as well as skyrmions in Bose-Einstein 
condensates (Ruostekoski and Anglin, 2001; Kasamatsu et at, 2005) and liquid crys- 
tals (Dicrking, 2003). Coherent polariton ensembles represent a promising new route 
to study the dynamics of defects and solitary waves (Amo et al., 2011). 

The study of multi-component fields is certainly among the most interesting new 
directions of research in this context. We have taken a first step in this direction by 
investigating the two-component Bose gas (Karl et al., 2013). The possibility of differ- 
ent non-thermal fixed points depending on inter- and intraspecies couplings opens a 
perspective on new types of experiments far from equilibrium (Kasamatsu and Tsub- 
ota, 2006; Nicklas et al., 2011). Relating (quasi-)topological field configurations known 
from equilibrium spinor Bose gases (Ueda, 2012) to transient scaling phenomena ob- 
served in correlation functions provides a great challenge for experiments and theory. 
Experimental studies of ultracold spin-1 and spin-2 Bose gases, including the detec- 
tion of spin domains, are far developed (Miesner et al., 1998; Schmaljohann et al., 
2004; Chang et at, 2004; Higbie et al, 2005; Sadler et al, 2006; Guzman et al, 2011). 
Multi-component fields are important far beyond ultracold atomic physics. For exam- 
ple, multi-component inflatons and their associated topological defects are discussed 
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in early-universe expansion scenarios (Rajantie, 2003; Berges, 2005; Allahverdi et al, 
2010). Exciting recent developments in the field of heavy-ion collisions, where the non- 
equilibrium regime of a quark-gluon plasma can be investigated by multi-component 
gauge field simulations (Arnold and Moore, 2006; Berges et al, 2009; Carrington and 
Rebhan, 2011; Fukushima and Gelis, 2012; Fukushima, 2011; Berges and Sexty, 2012; 
Berges et al., 20126) may exhibit a strong relevance of the type of interrelations dis- 
cussed above. Defect-induced non-thermal fixed points in this system are presently 
being explored. The rapid expansion of the quark-gluon plasma adds a completely 
new aspect to the dynamical description. Similar processes can be studied by releasing 
ultracold gases from their trapping potential, specifically in view of expanding turbu- 
lent clouds (Henn et al., 2009; Seman et al., 2011; Caracanhas et al, 2012; Weckesser, 
2012). 

Ultimately, non-thermal fixed points have to be included into a global picture 
of non-equilibrium dynamics of interacting many-body systems (Polkovnikov et al, 
2011; Gasenzer, 2009). The concept behind them points out a way towards universal 
phenomena far away from equilibrium, having many aspects in common with univer- 
sality and critical phenomena in thermal equilibrium. To understand their relations to 
non-thermal equilibrium states (Rigol et al., 2007; Eckstein and Kollar, 2008; Kollar 
and Eckstein, 2008; Rigol, 2009; Kronenwett and Gasenzer, 2011), generalised Gibbs 
states (Rigol et al, 2007; Eckstein and Kollar, 2008; Kollar and Eckstein, 2008), or 
prcthcrmalised states (Berges et al., 2004; Bonini and Wetterich, 1999; Aarts et al., 
2000; Barnett et al, 2011; Kitagawa et al, 2011; Gring et al, 2012; Kollar et al, 2011; 
Werner et at, 2012; Tsuji et al., 2012) are essential steps towards a unifying framework 
of complex dynamical many-body systems. To set up this framework, the development 
and extension of renormalisation-group techniques for far-from-equilibrium dynamics 
seems in order, and promising progress has been seen in the recent past (see Refs. given 
at the end of Sect. 2.2). 
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